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THE  NUMERICAL  SOLUTION  OF  COUPLED  NONLINEAR 
DIFFERENTIAL  EQUATIONS  OCCURRING  IN 
OCEANOGRAPHIC  THERMOSTAD  MODELING 


I.  INTRODUCTION 

Traces  of  temperature  as  a  function  of  depth  show  evidence  of  isothermal  masses  of  water  (ther- 
mostads)  embedded  in  an  otherwise  stratified  ocean.  These  thermostads  generally  reside  at  depths  of 
400-1000  m,  have  thicknesses  of  approximately  100  m,  and  diameters  no  greater  than  50  km  (see  also 
Broome,  Teague,  and  Hallock,  1979,  and  Rossby,  1980).  To  our  knowledge,  no  theoretical  work  has 
been  done  to  model  their  shape  and  velocity  fields;  although,  two  recent  Gulf  Stream  ring  papers  are  of 
related  interest. 

Csanady  (1979)  has  constructed  a  two-layer  model  of  warm  core  rings  employing  conservation  of 
potential  vorticity  and  radial  momentum.  Flierl  (1979)  then  expanded  this  work  by  including  the 
cyclostrophic  term  [(azimuthal  veiocity)2/radius]  which  arises  in  the  transformation  to  circular  coordi¬ 
nates.  These  principles  are  also  incorporated  in  this  work.  We  model  the  thermostad  as  an  oblate  oval 
solid  of  constant  density  water  (Fig.  1)  which  rotates  anticycionically  relative  to  a  stationary  ocean  hav¬ 
ing  a  linear  density  stratification.  We  then  seek  steady,  axisymetric  solutions  to  the  inviscid  equations 
of  motion  within  the  region  of  constant  density,  and  ignore  the  (2  effect  because  of  the  small  size  of  the 
thermostad.  The  equilibrium  shape  then  results  from  a  balance  between  the  pressure  gradient,  Coriolis 
and  cyclostrophic  forces. 


Manuscript  submitted  on  December  10.  1980. 


1 


MIGNEREY  AND  M1ED 


II.  EQUATIONS  OF  MOTION 


The  inviscid  Navier-Stokes  and  continuity  equations  in  a  rotating  cylindericai  coordinate  system, 
as  simplified  by  the  assumption  of  steady,  axisymetric  flow  inside  the  region  of  constant  density  <pn) 
are: 


—  +  /v  -  -L  /> 
r  Po 

(1) 

/«  -  0 

(2) 

P,  ~  ~PoS 

(3) 

v  -  v(r) 

(4) 

where  (u.  v)  are  the  velocities  in  the  (r.O)  directions  and  z  is  positive  upwards.  The  existence  of  Taylor 

d  9 

columns  can  be  shown  by  adding  the  quantities  —  (1)  and  —  (3).  Then 

dz  dr 

= - —  P„  =  0  which  proves  that  v  is  independent  of  z.  The  most  natural  wav  to  remove  the  degen- 

P  o  ' 

eracy  of  (2)  and  (4)  is  to  introduce  the  principle  of  conservation  of  potential  vorticity,  that  is 


at  +  / 

L 


constant 


where 


<u  ~  — 
r 


is  the  relative  vorticity,  /  the  Coriolis 


f  <»>-! 

dr  r  d9 


nnr'iiYiatur  iinrl  I 


(5) 

(6) 


To  obtain  a  pair  of  equations  for  the  velocity  and  column  height  it  is  necessary  to  express  the 
9  P 

radial  pressure  gradient  in  terms  of  the  column  height.  Taking  advantage  of  the  vertical  symmetry 
or 

of  the  oval  about  a  horizontal  plane,  we  introduce  the  half-height  \h  \  =  LI 2  so  that  ±h(r)  defines  the 
surface  of  the  thermostad.  Since  the  external  pressure  field  of  the  ocean  is  known,  it  is  possible  to  find 
the  pressure  at  h  ( r ) .  The  vertical  pressure  gradient  in  a  linearly  stratified  ocean  is 

|j--po*  (1-Sz). 

The  pressure  may  be  obtained  by  integrating  once: 

P{z)--Pr>g  (z-  1/2  8z2)  +  P0 
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where  Po  is  a  reference  pressure  at  r  -  0.  Likewise,  the  pressure  field  inside  the  thermostad  is 

P(r,z)  -  -pogz  +  c (r). 

Both  pressure  fields  must  match  over  the  entire  surface  of  the  thermostad 

P(h(r))  -  -  p0g  (h  -  1/2  bh2)  +  P0  -  -pogh  +  c(r) 

so  that 


Hence,  inside  the  thermostad 


c(r)  -  1/2  pog&h2  +  Pq. 


dP  ./^  «  6l/i2l 

—  -l/2Pos8  — 


(7) 


Substitution  of  (7)  into  (1)  produces 


v2  ,  ,  1  .  dh 2 

T  +  /v~  ST- 


(8) 


We  now  have  h  and  v  in  a  single  nonlinear  equation.  Another  equation  may  be  obtained  from  (5)  and 
(6): 


dr  r 


M 

Hn 


-  1 


(9) 


where  2 H0  is  the  length  to  which  a  Taylor  column  must  be  stretched  to  bring  the  relative  vorticity  to 
zero.*  We  look  for  a  family  of  solutions  corresponding  to  columns  which  have  contracted,  causing  a 
negative  relative  vorticity.  For  h  ~  Hq  there  is  very  little  contraction  and  only  a  small  amount  of  rela¬ 
tive  vorticity.  In  the  opposite  extreme  ( h  «  Hq ),  the  contraction  is  large,  as  is  the  relative  vorticity. 
The  boundary  conditions  for  (8,9)  are 


v(0)  -  0, 
h  (r0)  -  0. 


(10) 

(11) 


where  r0  is  the  radius  of  the  thermostad. 


Equations  (8,9)  with  boundary  conditions  (10,11)  are  nondimensionalized  by  letting 

r  -  Rr,  v-  F0v.  h  -  Hn  h  (12) 


*lt  will  emerge  in  Section  VII  that  all  observed  solutions  are  clustered  around  MO)  -  Hn.  In  effect,  then.  2 Hr,  is  very  nearly  the 
thermostad  height. 
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SO  that 


2  V2  ^  .  1  dh1 

e*  —  +  ev  -  — 


r 

d\ 

dr 


2  dr 

1/51-  l 


v(0)  -  0. 

6  (?o)  —  0 

v0  Ho 

where  e  =  —  is  the  Rossby  number  and  R  =*  —  lg&) 
RJ  J 


(13) 

(14) 

(15) 

(16) 


HpN 

I 


is  the  radius  of  deformation. 


Since  the  velocity  scale  is  chosen  as  the  velocity  at  the  edge, 

v  (f0)  -  1  and  K„-  «/V//0;  (17) 

then  the  nondimensionalization  (12)  takes  the  form: 


r-y  Ho  r.  v  =  eNH0  v.  h  -  H0h.  (18) 

The  system  of  Eqs.  (13-16)  constitute  an  eigenvalue  probelm  for  e.  Because  e  and  v  consistently  occur 
in  each  equation  only  as  a  product,  the  eigenvalue  problem  is  greatly  simplified.  After  introducing  the 
change  of  variables  v  =  ev  the  following  initial  value  problem  is  solved: 


vl  +  v-  1/2# 

r  dr 

(19) 

4  +  T-lAl-l 

(20) 

dr  r 

v(0)  -  0 

(21) 

0  <  6(0)  <  1. 

(22) 

These  equations  are  integrated  beginning  at  the  origin  and  continuing  until  h(r)  -  0.  This  determines 
r0  and  v(r0).  The  corresponding  eigenvalue  is  then  given  by  e  -  |v(r0)|.  Each  initial  value  6(0)  is 
thus  uniquely  related  to  a  particular  r0  and  e . 


III.  PROPERTIES  OF  THE  SOLUTION 

The  system  (19-22)  has  been  solved  in  detail  on  the  computer.  However  certain  properties  may 
be  deduced  using  analysis.  Since  both  (19)  and  (20)  are  even  in  6,  the  surface  of  the  thermostad  is 

indeed  symmetric  about  :  -  0.  Furthermore  at  r  *■  rn  we  expect  the  pressure  gradient  h  ~  to  be 


I 

! 

< 

I 
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bounded.  Since  )/(?g)  -  0  we  anticipate - -  «>;  the  solutions  support  this.  Below,  we  obtain  the 

dr 

shape  in  the  limiting  case  hi 0)  «  1,  and  derive  the  properties  of  the  thermostad  for  all  MO)  in  the 
vicinity  of  r  -  0. 

A.  Solution  for  h (0)  «  1 


An  analytic  solution  may  be  obtained  for  (19-22)  in  the  limit  h  —  0.  Then  (20)  becomes 


which  has  the  solution 


-4L  +  4~-  i. 

Jr  r 


v  -  -  r/2  + 

r 


The  initial  condition  (21)  forces  c.  =  0  and 


Using  this  expression  in  (19),  we  see  that 


v  (?)  -  -  r/2. 


dh 2  _  1  .2 

dr  2 


or 


h2(r)  -  -  r2  +  c2. 

At  the  outer  edge  of  the  thermostad,  r  **  r0  and  h  —  0;  thus  ci  =  rl  and 


h(r)  -  4  yjr o  -  r2. 


(23) 


(24) 


(25) 


The  surface  is  thus  given  by  Ah2  +  r2  —  /'g  which  shows  the  thermostad  to  be  an  oblate  ellipsoid  with 


/•0  ••  2A(0).  Furthermore,  differentiation  of  (25)  verifies  the  conjecture  — 

dr 


—  °° .  However 


(24)  shows  that  is  finite  at  ?0.  Although  we  have  only  shown  that  \h  — -|  <  »  for  h(Q)  «  1, 
dr  dr 


the  numerical  solutions  of  section  IV  support  the  veracity  of  this  assertion  for  all  M0). 
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B.  Power  Series  Solution  Near  r  —  0 

An  analytic  solution  near  the  origin  is  needed  for  a  successful  implementation  of  the  computer 
solution  in  the  following  section.  For  small  r  the  power  series  expansions 

v(r)  =  £  v„  ?n 

rr— 0 

and 

h(r)  =  £  hn  rn 

n— 0 

are  introduced  into  (19)  and  (20).  We  may  then  group  the  terms  according  to  orders  of  r. 


0  -7  : 
r 

Vo  “  0- 

v0=  0 

(26a,b) 

0(1): 

2v0v,  +  v0  -  M,. 

2v[  *  h0  -  1 

(27a.b) 

0(r): 

2v0v2  +  vj2  +  V[  =»  2 h0ft2  +  h\, 

<  ^ 

« 

fM 

i> 

(28a, b) 

The  importance  of  these  equations  lies  in  the  implications  to  be  drawn  from  the  0(1)  terms. 
.  -  /io —  1 

Since  v0  =*  0  and  h0  ^  0,  Eq.  (27a)  implies  h\  =  0  and  (27b)  yields  V|  =  — - —  while  (28b)  gives 
v;  —  0.  Thus,  for  small  r{  we  may  write 

v(? )  -  r  +  0  (r2)  (29) 

and 

h  (r )  —  h0  +  0(r2).  (30) 

We  employ  (29)  and  (30)  in  the  computer  solution  below. 

IV.  COMPUTER  SOLUTION 


Equations  (19)  and  (20)  are  a  pair  of  first  order  nonlinear  differential  equations  which  may  be 
rewritten  as: 


i£-ii  +  v  (3i) 

dr  r 
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=  02) 

Jr  r 

where 

-  _  1  *2 

g  =  y  h  ■ 

els 

Two  advantages  accrue  from  this  change  of  variables.  Near  the  edge  of  the  thermostad  is  well- 

c Jr 

behaved,  and  any  standard  scheme  for  solving  a  well-behaved  first  order  system  can  be  used  *  We  have 

written  a  subroutine  ca.iled  DIFF  (see  Appendix  A)  which  evaluates  (31)  and  (32)  in  a  form  usable  by 

the  integration  procedure  DIFSYS.  Two  minor  problems  arise.  The  integration  procedure  begins  at  the 

origin  where  the  1  /r  terms  in  (31)  and  (32)  are  infinite.  In  order  to  begin  the  integration,  we  employ 

»  ds  d  v 

the  power  series  solution  (28,29)  at  r  —  0  to  obtain  —  «*  0  and  — 

dr  o  dr 

The  other  problem  occurs  at  the  edge  of  the  thermostad.  Resolution  of  the  radius  r„  is  no  better 
than  the  step  size  Ar.  On  the  other  hand,  the  choice  of  a  small  A r  would  increase  the  computer  time 
prohibitively.  Fortunately  the  accuracy  of  DIFSYS  is  of  order  (Ar)2  (Bulirsch  and  Stoer  1966  b),  so 
that  a  small  Ar  is  not  needed  to  maintain  accuracy  of  h  and  v.  Hence,  a  large  Ar  is  used  until  the  edge 
is  detected  (i.e.,  when  g  becomes  negative).  At  this  point  the  subroutines  backtrack  to  the  last  known 
point  having  g  >  0.  A r  is  then  divided  by  ten  and  the  integration  proceeds  in  this  manner  until  Ar  is 
less  than  the  convergence  tolerance.  The  evaluation  of  r0  is  then  considered  to  be  finished. 

Although  the  overall  accuracy  of  the  method  is  0[(Ar)2],  a  few  sample  calculations  were  made  to 
ensure  the  proper  choice  of  Ar.  The  results  summarized  in  Table  1  indicate  a  value  of  A r  =  1CT:  is 

,  Table  1 

Values  for  h(r)  and  v(r)  at  r  —  1 
for  various  values  of  Ar. 


1? 

Ad) 

v(l) 

#  of  computations 

ir* 

.22369 

-.30885 

10 

ir2 

.22357 

-.30886 

100 

lO-3 

.22256 

-.30897 

1000 

"tn  this  paper  we  use  the  routine  DIFSYS  which  is  available  through  the  NRL  Research  Computation  Division  It  is  an  extrapo¬ 
lation  method  based  upon  the  works  of  Bulirsch  and  Stoer  (1966  a.b) 
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appropriate  for  three  place  accuracy  in  h  and  v.  The  smaller  value  requires  too  many  calculations  and 
the  larger  value  produces  too  few  data  points.  Accordingly  the  convergence  tolerance  for  1?  was 
chosen  to  be  10_J. 

v.  numerical  results 


The  functions  hi,?),  v  (?)  and  the  radius  ?0  were  computed  for  values  of  h  (0)  =  .05,  .10, 

.90,  .95,  1.00,  1.05.  The  results  are  summarized  in  Figs.  2  and  3.  Figure  2  shows  hi?)  for  several 
values  of  hi 0);  these  curves  correspond  to  different  values  of  «.  Closed  surfaces  exist  when  MO)  <  1; 
note  that  no  such  surface  exists  for  the  curves  M0)  =  1.00,  1.05.  The  associated  velocity  field  (Fig.  3) 
is  negative  and  monotonically  increasing  in  magnitude  more  rapidly  than  a  linear  function  of  ?. 
Although  the  radial  velocity  gradient  is  larger  for  small  M0),  the  large  size  of  ?o  for  large  M0)  pro¬ 
duces  a  larger  peak  velocity. 


For  small  values  of  e 


NH0 


the  computer  solutions  approach  the  asymptotic  solutions  (23) 


*  -  Vh  (0) 

and  (25).  In  particular,  the  aspect  ratio  a  =  hi 0)/?0=  -  approaches  0.5,  as  shown  in  Fig.  4. 

Jro 

For  larger  values  of  €  the  aspect  ratio  is  smaller,  which  indicates  the  thermostad  is  more  oblate  when 
the  fluid  has  less  relative  vorticity.  This  raises  the  question  as  to  how  much  the  function  hi?)  deviates 
from  the  ellipsoidal  limit  («  ~  0)  when  e,  M0)  and  ?0  assume  larger  values.  To  answer  this,  the  vari¬ 
ables  were  normalized  as  follows: 


r*  -  ?/?q,  h'/hiO),  v*=*  v/v(f0). 

As  a  result  of  this  normalization,  all  thermostads  transform  to  shapes  with  unit  radius  and  unit  height. 
Some  normalized  thermostads  are  shown  in  Figs.  5  and  6  for  two  different  values  of  « .  For  the  larger 
value  of  €  the  surface  of  the  thermostad  is  more  square  than  elliptical,  and  the  velocity  is  markedly 
nonlinear. 


As  indicated  above,  two  pathological  cases  were  run.  For  M0)  »  1.00  the  Taylor  columns  exist  at 
their  initial  (pre-collapse)  height;  thus  there  is  no  relative  vorticity.  Hence  the  velocity  is  zero  every- 
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where  and  the  bounding  thermostad  surface  extends  to  infinity.  When  M 0)  >  1,  e.g.,  /;(0)  =  1.05. 
the  Taylor  columns  are  stretched  producing  positive  relative  vorticity  and  corresponding  positive 
(cyclonic)  velocities.  The  surface  is  again  unbounded  as  it  rises  away  from  the  origin. 

VI.  A  METHOD  FOR  OBTAINING  THE  ROSSBY  NUMBER  AND  HEIGHT  SCALE 

The  shape  and  velocity  field  of  a  thermostad  depend  on  the  value  of  e.  However  e  is  not  easy  to 
determine  for  a  dimensional  thermostad  even  if  v(r0)  is  known.  This  is  because  the  scale  height  //n 
appears  in  all  of  the  nondimensionalizations  (18),  and  Hn  is  not  directly  measureable  for  a  thermostad 
in  the  ocean.  Two  measurements  which  are  less  difficult  to  make  are  /t(0)  and  r0.  Since  each  of  these 
quantities  depends  on  the  Rossby  number,  the  aspect  ratio  a  (-  MO)//*,,)  is  also  a  function  of  €, 
a  =■  a(e).  This  functional  dependence  is  shown  in  Fig.  4.  From  the  nondimensionalization  (18),  we 

may  express  aasn*-1-  h_i0)_  gjnce  t^e  ^gp^  /,(Q)  an(j  radjus  r0  can  be  estimated  from  the  data,  € 
J  ro 

may  be  found  for  and  observed  thermostad.  Then  the  height  scale  is  given  by  ff0  -  h( 0)/M0;e) 
where  MO;*)  is  shown  in  Fig  7a;  the  results  for  ra  as  a  function  of  e  are  shown  in  Fig.  7b.  With  *  and 
Ht)  thus  obtained.  Figs.  2  and  3  and  the  scaling  (18)  may  be  used  to  estimate  the  shape  and  velocity 
field  in  dimensional  form. 

One  further  problem  occurs  when  the  data  is  reduced.  The  radius  r0  is  not  precisely  known;  only 
an  upper  bound  can  be  established  by  examining  contiguous  XBT's.  This  leads  to  inequalities  for  a.  e. 
Ho  and  E0.  In  trying  to  establish  the  inequality  for  F0,  Eq.  (17)  is  needed.  Unfortunately  the  inequali¬ 
ties  are  such  that  *  is  less  than  some  bound  and  H0  greater  than  another  bound,  leaving  the  product 
indeterminate.  Likewise  MO;*)  =  MO)///0  implies  M 0;e)  is  less  than  a  bound.  To  avoid  this  prob- 

*  vo 

lem,  a  velocity/height  ratio  is  introduced.  /3  =  — -  —  —  —  (33).  This  ratio,  which  is  a  function 

MO;*)  <VM0) 

of  *.  is  shown  in  Fig.  8.  The  slope  of  this  curve  in  conjunction  with  the  inequality  for  *  implies  (3  is 
greater  than  a  bound.  Likewise  the  equation  E0  =■  /3,VM 0)  means  the  peak  velocity  can  only  be  esit- 
mated  as  exceeding  some  minimum  value. 
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VII.  APPLICATION  TO  OBSERVED  THERMOSTADS 

Traces  of  temperature  as  a  function  of  depth,  obtained  by  several  XBT’s,  show  evidence  of  ther- 
mostads  (Figs.  9a-g).  Aside  from  Figs.  9d  and  9f,  however,  there  is  only  a  single  trace  through  any  of 
the  suspected  thermostads.  In  no  case  is  a  velocity  measurement  available;  consequently,  an  estimate 

of  the  aspect  ratio  a  —  ■—  -  ^  ,  peak  velocity  v(r0)  and  height  scale  /i(0)  is  all  that  can  be  obtained. 
J  | 

To  accomplish  this,  the  Brunt-Vaisala  frequency  must  be  evaluated.  Based  on  historical  temperature- 
salinity  data  corresponding  to  the  location  or  the  XBT,  the  salinity  and  rrT  are  computed  as  a  function 
of  depth.  The  Brunt-Vaisala  frequency  is  then  obtained  using  the  central  difference  formula 

N*  =  —  -  at  depth  z  —  z0-  Since  differentiation  produces  varying  results  due  to  fluctuations 

Po  [  -V  ~  J- 

in  the  data,  we  formed  two  averages  for  N:  one  actually  within  the  thermostad.  and  the  other  over  the 
surrounding  neighborhood. 

The  local  value  for  /was  found  using  the  known  latitude  of  the  XBT.  The  vertical  location  of  the 
thermostad  was  obtained  by  drawing  a  straight  line  through  the  thermostad  and  subjectively  determin¬ 
ing  where  the  temperature  gradient  changes  slope.  Then  the  thickness  was  measured  from  the  plot  and 
halved  to  form  h(0).  Because  of  the  sparseness  of  the  data,  only  an  upper  bound  can  be  set  for  r0  based 
on  the  distance  between  contiguous  XBT’s.  Using  the  data  in  conjunction  with  the  theory  allows 
bounds  to  be  put  on  « ,  a ,  H0  and  V0. 

As  summarized  in  Table  2,  the  observed  thermostads  are  oblate,  having  an  aspect  ratio  in  the 
range  .119  <  a  ^  .182.  This  low  value  of  a  is  due  to  the  thinness  of  the  observed  mass  of  water  in 
conjunction  with  the  large  upper  bound  on  the  diameter.  Although  the  thermostads  are  thin,  the 
heights  of  their  Taylor  columns  are  not  much  less  than  the  scale  height.  Indeed,  h  (0)  which  is  a  meas¬ 
ure  of  this  contraction  is  in  the  range  .992  to  .999.  A  significant  result  of  this  work  is  H0  approximately 
equals  half  the  thermostad  height,  with  less  than  1%  error. 
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Table  2 

A  Comparison  of  the  Observed  Parameters  and  Derived 
Inequalities  for  Several  Thermostads 


Figure 

4a 

9b 

9c 

9d 

9e 

or 

9g 

'0 

(km) 

26 

22 

26 

I- 

27 

28 

32 

- - - j 

25 

hi 0)  im>  * 

i  *» 

65 

65 

85 

100 

47 

52 

v 

T  outside 

rad 

s 

- 

V 

o 

VS 

oo 

r*~> 

.366  10'-’ 

.403  10'2 

.373  ■  10--’ 

.407  I0~2 

430  I0~2 

.379  I0-2 

^  inside 

rad 

* 

- 

.310  •  10~2 

190  10~2 

124  10' 2 

.093  ■  10‘2 

.244  ■  I0-2 

324  ■  10*2 

.317  •  nr2 

t 

rad 

$ 

,862  10-4 

.885  10"* 

.834  10"* 

.794  10"* 

.800  -  10"* 

.899  •  I0-* 

.876  lO*4 

f  « 

> 

.167 

122 

.121 

.148 

182 

070 

090 

€ 

$ 

553 

561 

.561 

557 

551 

566 

565 

a 

559 

.561 

561 

.560 

559 

566 

;65  ; 

0 

< 

6.04 

8.232 

8.232 

6.930 

5.472 

11  294 

10.680 

*0 

< 

.992 

999 

999 

995 

.987 

9999 

9999 

Hn  (ml 

r 

101.3 

47 

52 

LJji 

(m/s)  ^ 

-23g 

114 

Ml  | 

The  observed  thermostads  appear  to  have  a  weak  relative  vorticity.  The  Rossby  number  which  is 
obtained  from  Fig.  4  once  a  is  known  lies  in  the  interval  .551  to  .561  and  is  large  enough  for  significant 
departures  from  geostrophic  flow  to  occur.  An  estimate  of  the  peak  velocity  may  be  obtained  using 
(33).  Based  on  the  bounds  for «,  Fig.  8  is  used  to  obtain  /3,  which  lies  in  the  range  .559  to  .562.  The 
associated  peak  velocities  are  then  .132  (m/s)  to  .228  (m/s).  Although  these  velocities  are  large,  we 
believe  they  are  reasonable  in  view  of  the  theory’s  neglect  of  a  Brunt-Vaisala  frequency  inside  the  ther- 
mostad  (see  Table  2).  An  inclusion  of  non-constant  N (: )  would  lead  to  smaller  velocities. 


Clearly,  a  measurement  of  the  velocity  field  inside  a  thermostad  would  be  a  more  stringent  test  of 
the  theory.  Unfortunately  such  measurements  may  be  difficult  to  make.  A  more  fruitful  approach 
would  entail  dropping  a  sufficient  number  of  XBT's  to  enable  a  comparison  of  the  predicted  and 
observed  shapes,  along  with  a  better  determination  of  r0. 
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